%% 清空变量环境
close all;
clc,clear;
warning off;

%%  导入数据
res=readmatrix("sample_random_forest_data.xlsx");
res1=res(991:end,1:end-1);%新预测数据
res=res(1:990,:);%样本数据
%% 数据归一化(去除量纲的影响)
X=res(:,1:end-1);
Y=res(:,end);
[x,psin]=mapminmax(X',0,1);

%% 划分数据集与测试集
num=size(res,1);% 总样本数
k=1;
if k==0
    statr=1:num;%不打乱样本
else
    statr=randperm(num);%打乱
end
r=0.8;%训练集比例
trainnum=floor(num*r);
xtrain=x(:,statr(1:trainnum));
xtest=x(:,statr(trainnum+1:end));
Ytrain=Y(statr(1:trainnum));
Ytest=Y(statr(trainnum+1:end));
%%训练模型
trees=100;  %决策树个数
leaf=5;    %叶子数，过小会过拟合
wucha='on'; %打开误差图
importance='on';    %计算特征的重要性
net=TreeBagger(trees,xtrain',Ytrain','OOBPredictorImportance',importance,'Method','classification','OOBPrediction',wucha,'MinLeafSize',1,'MaxNumSplits',50,'minleaf',leaf);
import=net.OOBPermutedPredictorDeltaError; %重要性

%% 预测
[re1,scores1]=predict(net,xtrain');
[re2,scores2]=predict(net,xtest');
%格式转化
pre1=str2num(cell2mat(re1));
pre2=str2num(cell2mat(re2));

%% 数据排序
[ytrain,idx1]=sort(Ytrain);
[ytest,idx2]=sort(Ytest);

Pre1=pre1(idx1);
Pre2=pre2(idx2);

%% 误差
error1=sum((pre1==Ytrain))/trainnum*100;
error2=sum((pre2==Ytest))/(num-trainnum)*100;

%% 混淆矩阵
figure
cm=confusionchart(ytrain,Pre1);
cm.Title='训练集混淆矩阵';
cm.ColumnSummary='column-normalized';
cm.RowSummary='row-normalized';

figure
cm=confusionchart(ytest,Pre2);
cm.Title='测试集混淆矩阵';
cm.ColumnSummary='column-normalized';
cm.RowSummary='row-normalized';
%% 新数据预测
newdata=res1;
[newy,scores]=newpre(newdata,net);
%保存数据
% xlswrite()

%% 画图
%误差图
figure
plot(1:trainnum,ytrain,'r-^',1:trainnum,Pre1,'b-*','LineWidth',1)
legend('真实值','预测值')
xlabel('样本点')
ylabel('预测值')
title(['训练集预测结果对比','准确率:',num2str(error1)])

figure
plot(1:num-trainnum,ytest,'r-^',1:num-trainnum,Pre2,'b-*','LineWidth',1)
legend('真实值','预测值')
xlabel('样本点')
ylabel('预测值')
title(['测试集预测结果对比','准确率:',num2str(error2)])


%误差曲线(带外错误率曲线)
figure
plot(1:trees,oobError(net),'r--o','LineWidth',1)
legend('误差迭代曲线')
xlabel('决策树(迭代次数)')
ylabel('误差')
xlim([1,trees])
title('误差迭代曲线')

%特征重要性图
figure
bar(import,'green')
yticks([])
xlabel('特征')
ylabel('重要性')

%ROC曲线
% 获取类别名称
classNames = net.ClassNames;

% 初始化图形
figure;
hold on;
grid on;

% 为每个类别计算和绘制 ROC 曲线
for i = 1:numel(classNames)
    positiveClass = classNames{i};
    positiveClassIdx = strcmp(net.ClassNames, positiveClass);
    scoresPositive = scores1(:, positiveClassIdx);
    
    % 计算 ROC 曲线数据
    [Xroc, Yroc, ~, AUC] = perfcurve(Ytrain, scoresPositive, positiveClass);
    
    % 绘制 ROC 曲线
    plot(Xroc, Yroc, 'LineWidth', 2, 'DisplayName', [positiveClass ' (AUC = ' num2str(AUC) ')']);
end

% 设置图形标签和标题
xlabel('假正利率(FPR)')
ylabel('真正利率(TPR)')
title(['ROC Curve (AUC = ' num2str(AUC) ')'])
legend show;
hold off;

%% 预测函数
function [y,scores]=newpre(newdata,net)
    x=newdata;
    %预测
    [ty,scores]=predict(net,x);
    %格式转化
    y=str2num(cell2mat(ty));
end

